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Abstract 

An Ornstein-Zernike approximation for the two-body correlation function em- 
bodying thermodynamic consistency is apphed to a system of classical Heisenberg 
spins on a three-dimensional lattice. The consistency condition determined in a 
previous work is supplemented by introducing a simplified expression for the mean- 
square fluctuations of the spin on each lattice site. The thermodynamics and the 
correlations obtained by this closure are then compared with approximants based on 
extrapolation of series expansions and with Monte Carlo simulations. The compar- 
ison reveals that many properties of the model, including the critical temperature, 
are very well reproduced by this simple version of the theory, but that it shows sub- 
stantial quantitative error in the critical region, both above the critical temperature 
and with respect to its rendering of the spontaneous magnetization curve. A less 
simple but conceptually more satisfactory version of the SCOZA is then developed, 
but not solved, in which the effects of transverse correlations on the longitudinal 
susceptibility is included, yielding a more complete and accurate description of the 
spin-wave properties of the model. 

1 Introduction 

Several approaches in the context of liquid-state theory hinge on the Ornstein-Zernike 
(OZ) equation that relates the two-body distribution function g{r) and the direct correla- 
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tion function c(r) . Once an approximate expression for c(r) is inserted into this equation, 
a closed theory for the thermodynamics and the correlations is obtained. 

While the dependence of c(r) on the thermodynamic variables is fed directly into the 
expression of c(r) itself in the commonly used closures to the OZ equation, another method 
for obtaining a closed theory amounts to requiring that c(r) depends on the density and 
on the temperature of the system in such a way that consistency between different paths 
to the thermodynamics is obtained. Recent results on simple lattice gases |l| and model 
fluids have shown that such a self- consistent Ornstein-Zernike approximation (SCOZA) 
is able to yield very good results for both the overall thermodynamics and correlations 
and for the critical point and the coexistence curve. 

The idea of using thermodynamic consistency as a constraint to close the OZ equation 
is relevant not only to simple fluids. In particular, the theory has already been formu- 
lated also for the vector model 0, i.e., a system of classical D-dimensional spins with 
nearest-neighbor interaction. In the present paper the D-vector model SCOZA is fur- 
ther investigated and the simplest version of the SCOZA that we have developed is solved 
numerically in the case of three-dimensional spins (classical Heisenberg model) on a three- 
dimensional lattice, and the results are compared with Monte Carlo simulations [^, ^, § 



and Fade approximants based on high-temperature expansions [0, |, [T^, |Tl| . An anal- 
ysis of what is lacking in this simple application to the D-vector model then guides us in 
formulating a more sophisticated version of the theory that takes into account the influ- 
ence of transverse correlations on the longitudinal susceptibility and as a result, yields a 
more complete and accurate description of the spin-wave properties of the model. 

Besides allowing us to test the SCOZA for a system with a vector order parameter, 
we believe that the investigation set forth here is also interesting in its own sake. While 
the main properties of the D-vector model have been understood for a long time and 
its critical constants have been determined with great accuracy in several recent papers 



both by finite-size scaling techniques Q and by extrapolation of series expansions |]TT 
there is nevertheless a paucity of treatments that yield a reasonable description of the 
model for D = 3 (the classical Heisenberg model) over the whole phase plane, especially 
compared to the Ising or D = 1 case. The D-vector model has recently been studied 
by cumulant expansion methods [|T^ |TB|], which can be regarded as a systematic way 
to improve the mean-field approximation. However, results were determined only in zero 
magnetic field, and do not include the correlation functions. Moreover, the predictions for 
the critical temperature obtained within the currently adopted truncation of the cumulant 
expansion are in only modest agreement with the expected values |jl3| . It seems therefore 
interesting to develop a version of SCOZA for this model, particularly for D = 2 and 
D = 3, that promises to yield high overall accuracy, especially since our approach gives 
both the thermodynamics and the correlations of the system. 

A difference that arises in the present implementation of the SCOZA with respect to 
the Ising case is that the OZ approximation (i.e., the assumption that the direct correla- 
tion function c(r) has the same range as the interaction between the spins) and the short- 
range constraint on the two-body correlation function which stems from the fixed length 
of the spins are not enough to get a closed theory by means of a single thermodynamic 
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consistency condition. We need some ansatz about the mean-square fluctuations of tlie 
transverse and longitudinal spin components. The initial choice we investigate here is to 
regard these fluctuations of the spin on each lattice site as a function of the magnetization 
alone, thereby using at every temperature the expression valid in the high-temperature 
limit. We find that within this approximation the error in the critical temperature is 
at most about 0.5%, and the compressibility, the correlation length, and the correlation 
functions are nicely reproduced over a wide range of thermodynamic states. Nevertheless, 
the global accuracy of this version of the theory is inferior to that of the Ising or simple 
fluid case. In particular, for both the Ising and the Heisenberg model the SCOZA gives 
spherical critical exponents above the critical temperature Tc, but in the Ising case 
the asymptotic critical regime is detectable only very close to the critical point and the 
thermodynamics is well reproduced up to reduced temperatures t = (T — Tc)/Tc ~ 10~^. 
For the Heisenberg model, on the other hand, the deviations from the correct critical 
behavior appear already for t ~ 10^^. At smaller values of t, the theory readily reaches 
the asymptotic spherical regime. Moreover, while in the Ising-model SCOZA the coex- 
istence curve displays a non-classical exponent (3 = 7/20 |14|, the Heisenberg coexis- 
tence curve predicted by the SCOZA instead is associated with /5 = 1/2, and is overall 
considerably narrower than the expected one as given by MC simulations and Pade 
approximant techniques 0]. (We find however that it displays the correct behavior in 
the low-temperature limit). In order to probe the sensitiveness of the theory to the ap- 
proximate ansatz adopted for the mean-square spin fluctuations, we also implemented a 
different approximation, which was originally proposed in Ref. The results are almost 
identical in the two cases. In response to the lack of high quantitative accuracy shown 
by these simple versions of the theory we then develop (but do not solve here) an im- 
proved version that incorporates a more complete and accurate description of spin- wave 
correlations. 

The paper is organized as follows: in Sec. 2 we review the theory already developed 
in Ref. , introduce an approximation necessary to obtain a closed equation, and cast it 
into a form amenable to numerical solution, which we obtain in Sec. 3. There we show our 
results for the critical point, the thermodynamics, the correlations, and the phase diagram, 
and compare them to predictions from MC simulations and Pade approximants. We give 
our improved version of the theory in Sec. 4, in which we include the effects of transverse 
correlations on the longitudinal susceptibility. Our findings are briefly summarized in 
Sec. 5. 



2 Theory 

We consider a system of D-dimensional vectors Sr of unit length on a three-dimensional 
lattice, r being the position of a generic lattice site, interacting via a nearest-neighbor, 
ferromagnetic coupling. Here we are interested in the classical analog of the Heisenberg 
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model, and therefore will deal with the case D = 3. The hamiltonian of the system is 



n 



-J Sr 

(r,r') 



■Sr' 



HE 



(1) 



where the sum over r and r' extends over all nearest-neighbor couples, H is the external 
magnetic field, and J > is the coupling constant, which will be taken equal to unity in 
the following. The SCOZA for this model was originally introduced in Ref. 0. We begin 
by reviewing it in some detail. 

The reduced longitudinal and transverse isothermal susceptibilities xy, X± defined 

by: 
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where m is the magnetization per spin and (3 = l/k^T, T being the absolute temperature, 
and fee the Boltzmann's constant. (We remark that the definitions of X± used here 
differ from those in Ref. by a factor (3). Here the symbols || and _L refer respectively 
to the direction along the magnetic field H, and to any of the D — 1 directions of the 
coordinate axes lying in the plane perpendicular to H. Because of the spherical symmetry 
of the system, we find 
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where a = a{m, f3) is the Helmholtz free energy per spin. Let us note that Eqs. 
imply 

1 d ( m\ 
XII dm \x±) ' 
The internal energy per spin u is given by 



u 



We now introduce the two-spin correlation function matrix r(r) with elements 




(6) 



(7) 



r,,(r) = (44) - (4) (4) 



(8) 
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where the indexes i, j refer to two generic components of the vectors Sr, Sq. Isotropy 
imphes that T can be put into diagonal form by choosing the direction of the magnetic 
field H as one of the coordinate axis. In particular, we have 

rii(r) = (44)-m^ (9) 
rx(r) = (s^s^). (10) 

We recall that both the internal energy u and the susceptibilities xy, X± can be expressed 
by the correlation functions rii, F^. Specifically, we have for the internal energy 



u = — 
2 



+ T\\{X) + dT^{\) . (11) 



Here A is a vector joining a given lattice site with any of its nearest neighbors, q is the 
number of nearest neighbors, and d = D — 1. (In Ref. the factor q is not used. This 
corresponds to setting Jq = 1). The susceptibilities are given by 

XII = f||(/. = 0), (12) 
X± = f±(A; = 0), (13) 

where r||(k), r_L(k) are the Fourier transforms of r||(r), r_|_(r). Since the length of each 
spin is equal to unity, we have for every site r 

{{sir)+d{{s^f) = i. (14) 

If we set 

=r||(r = 0), (15) 

R± = r^{r = 0), (16) 

we find from Eqs. (|^), ( [T0| ) that the condition ( p!^ becomes 

Rll+dR± = l-m^ . (17) 

Let us now consider the direct correlation function c(r) whose Fourier transform c(k) is 
related to r(k) by 

where we have understood the subscripts ||, _L in both r(k) and c(k). In order to simplify 
the notation, we will adopt this same convention in the remaining of the Section whenever 
a relation applies to both the longitudinal and the transverse components of the correla- 



tion functions. Eq. (|18|) corresponds to the OZ equation commonly used in liquid-state 



theory [jT5|. We note that, since the matrix r(r) is diagonal, Eq. ([T8|) is a purely scalar 



relation. Here we use the ansatz also known as OZ approximation, namely we assume 
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that c(r) has the same range as the interaction between the spins. For the model in study 
we have then 

Co r = , 

c(r) = < ci r = A, (19) 
otherwise . 



Eq. ([T9|) is readily written in Fourier space: 

c(k) = Co + gci7(k) , 

where the nearest-neighbor Fourier transform 7(k) is given by 

1 



7(k) 



- E e 

^ A=n.n. 



(20) 



(21) 



the sum running over all the nearest neighbors of a given lattice site. We can now use the 
OZ equation (|18|) to express r(r) via c(k). Let us introduce the quantity z given by 



qRci 
. - Rcn 



and the lattice Green's function 

P{z) 



1 



B(27r)3 1 - ;Z7(k) 



(22) 
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where B is the first Brillouin zone of the lattice. The function P{z) has been widely 
studied in the literature |T^, [l^. For the simple cubic (SC), body-centered cubic (BCC), 
and face-centered cubic (FCC) lattices considered here it can be evaluated in terms of 
elliptic integrals [1^. /^From Eqs. ([T5|)-(|23|) we readily find 
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Cl 
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zP{z) 
qR 
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1 - P{^) 
zP{z) 



(24) 
(25) 
(26) 



The susceptibilities and the internal energy can then be expressed in terms of z and R. 
From Eqs. (0), (|ll), (0), (H), (H) we get 



1 

X 



while Eqs. (0) and (|6D yield 
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-m^ + R\\ 



R 



^||P(;.||) +^''^ZxP(^x) 



(27) 



(28) 
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It appears then that within the OZ approximation the thermodynamics and the corre- 
lations are determined by the quantities z\\, z±, R\\, R±. These quantities are not all 
independent: in fact, they must be chosen in such a way that the conditions (|) and 
(P!7D are satisfied. As noted above, the former is a consequence of the isotropy of the 
system, while the latter stems from the requirement of each spin having unit length. The 
other constraint we want to be satisfied is that of consistency between different routes 
to the thermodynamics which is peculiar to the SCOZA. According to the energy route, 
the thermodynamics is determined from the internal energy of the system as given by 
Eq. (11). On the other hand, knowledge of the correlations may be used to obtain the 
thermodynamics also via Eq. (|T^) or (|13]) for the susceptibilities. This is the analog of 



the compressibility route in fluid systems. While in a hypotetical exact treatment both 
routes are obviously bound to yield the the same result, this is general no longer true as 
soon as approximations are introduced. As is readily found from Eqs. (|p, (^, (|^), the 
requirement that the susceptibilities X± ^-^id the internal energy u come from a unique 
Helmholtz free energy is embodied in the relations 



(29) 
(30) 
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We note that these expression are not independent as they are related by Eq. (H). Sub- 
stitution of Eqs. (|2^), (|28|) into Eq. (BO) yields the thermodynamic consistency condition 



1 A 



1 + — — 

2m 8m 



Rn 



1 - P{Z\\ 

^ll^(^ll) 



+ dR_ 



1 - Piz±) 
z±P{z±) 



(31) 



and there will be a corresponding equation obtained from Eq. (p9|) . Eqs. @, ([TTD , and 
( |3TD are a set of three independent conditions for z\\, z±, R\\, R±. As we have four unknown 
quantities, this is clearly not enough to obtain a closed set of equations. We then need 
one more equation, which will be introduced as a further approximate ansatz besides the 
OZ approximation (|19[). This is different from what is found in the Ising model or D = 1 
SCOZA, in which the OZ ansatz is the sole approximation introduced in the theory. The 
closure relation suggested in Ref. was 

4 = (32) 

i.e., the nearest-neighbor contributions to the direct correlation function in the longitudi- 
nal and in the transverse directions are assumed to coincide. Eq. (|25| ) then gives 



Ri 



R 



z^Pjz^) 



(33) 



On the other hand, Eq. (^7\ 
susceptibility corresponds respectively to z\\ 



shows that a divergence of the longitudinal or transverse 

T§, in the D- 



1 or z_L = 1 . As is well known 



vector model on a three-dimensional lattice with D > 2 both susceptibilities are expected 
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to diverge on the spontaneous magnetization curve, //this is actually the case also in the 
present approximation, then Eq. ( |55D would imply that on the spontaneous magnetization 
curve the longitudinal and transverse mean-square deviations of the spin coincide. This 
seems quite unphysical, especially at low temperature, when one expects the behavior of 
the system to be governed by small transverse oscillations around the direction of the 



average magnetization m |T9|, thereby yielding R\\ <^ i?^. As a consequence, although 
we did implement the approximation of Eq. (^) on order to probe the sensitiveness of 
the SCOZA to the closure adopted (see Sec. 3), we decided to resort to a different one 
for most of the computations performed here. Specifically, we simply disregarded the 
dependence of R\\, R± on the temperature of the system, i.e. we took for these quantities 
the same espressions one finds in the high-temperature limit. Since in this limit the 
spins on different sites are uncorrelated, it is found from Eqs. (0), (0) that the present 
approximation amounts to assuming, both for the parallel and the transverse components 
of R, 

R = x{m,P = 0). (34) 

The reduced susceptibilities in the high-temperature limit are just those of the non- 
interacting spin system. They are then straightforwardly determined from the Langevin 
function: ^ ^ 

= tanh{(3H) ~ Jh ' ^^^^ 



where the subscript "id" refers to the non-interacting system. Eq. ( |5BD is then inverted to 
express R\\) R± as a function of m. We observe that in the approximation of Eq. ( ^l)) the 
condition (|T7p is trivially satisfied at every temperature. Eqs. (^ and (|S^ ) finally allow 
one to obtain a closed partial differential equation (PDE) for the quantity z± from the 
thermodynamic consistency condition (|3T|). We now perform some manipulations on this 
equation in order to cast it into a form suitable for numerical integration. Let us note 
that Eqs. (|25|) , (|27|) imply < z < 1. In this interval the function 

y={l- z)P{z) (36) 

is invertible. We then introduce its inverse function z = F{y), < y < 1, and we set 

G(y) = '-y~^^y^ - l-^(^) (37) 
yF{y) - zP{z) ' ^^^^ 



R\\G{y\{) + dR^G{y^) . (3 



so that Eq. (plf ) becomes 



1 dy± _i : ^ ^ 



qR± dp 2m dm 

We now use Eq. (^ along with Eqs. (p?!), ( p6D to express in terms of y±: 
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In the present approximation in which R± coincide with the longitudinal and trans- 
verse susceptibilities of the non-interacting system, it is readily found using Eq. (^) that 
Eq. (|39|) becomes 

y\\ =y± + m— . (40) 



Moreover, we find 



dy\\ / m ^ 
dm \ R\ 



R I dm 



dm R I dm^ 



(41) 



where we have indicated with R'\\ the first derivative of R\\{m). Eqs. (^0]), (|T) are then 
substituted into Eq. (|38|) to yield 



2 dy, 



^P'^y\\)-Q^ + {^^ [^R\\G\y\{) + dR^G'iy^)] + mi?||i?'||G'(2/||)} 



g d(3 



+ 



'R!\\G{y\\)+dR!^G{yi_)-2m 



1 



m 9m 
(42) 



Since R\\, R± are known functions of m and y\\ is given by Eq. (|40|) , the above equation 
is a closed PDE for as a function of /? and m. Let us observe that Eq. (^) and 
Eq. ( ^71) imply that both ?/|| and y± are even functions of m. As a consequence, the 
term m'^dy±/dm which appears in the r.h.s. of Eq. ( |^ ) does not imply any singularity 
at m = 0, since in this limit it just gives d'^y±/ dm^ . The PDE (^2]) is of course non- 
linear. In particular, both y± and dy±/dm (via Eq. (|40|) ) appear as arguments of the 
functions G, G' . However, the dependence on d'^y^jdm? is linear. When solving Eq. ( ^2]) 
numerically, one can then take advantage of this structure by resorting to some of the 
finite-difference algorithms especially devised for gfiasz-Zmear diffusive PDE's |2^. While 
these algorithms share with fully implicit methods the unconditional numerical stability, 
they do not require the solution of a large system of coupled non-linear equations at each 
integration step and are therefore easier to implement. 

In order to integrate Eq. ( ^2]) we also need a set of initial and boundary conditions. 
The initial condition is given at /5 = 0, i.e. in the high-temperature limit. In this limit all 
the correlations vanish except onsite so that, as we pointed out above, x± and R^ coincide. 
Eq. (|27|) implies then = 0, which in turn gives via Eq. (RBI) the initial condition 



yi_{m,i5 = ^) 



for every m. 



(43) 



The magnetization per spin m varies in the interval — 1 < m < 1. However, since yx_ is an 
even function of m, we can restrict the integration of Eq. (^) to the interval < m < 1 
provided we require y^ to be indeed symmetric at m = 0: 



dyi. 
dm 



(m = 0,/?) = for every /?. 



(44) 



As m approaches its limiting value m = 1, the longitudinal and transverse susceptibilities 
X\\ and x± vanish respectively as (1 — m)^ and 1 — m. In this limit the reciprocal suscep- 
tibilities diverge like those of the non-interacting system, as the excess contribution stays 
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finite. We have then again ~ V-^J- hence from Eqs. (pTf) , (p6D 



l/_L(m = l,/5) = 1 for every (3. (45) 

We note that for the transverse susceptibihty to be positive, y_[_ has to be positive as well. 
The numerical solution of Eq. ( ^21) shows that this condition is always satisfied if the 
inverse temperature [3 is smaller than a certain value (3c such that y±(m = 0,/5c) = 0. At 
13 = (3c both the transverse and (from Eq. (|^)) the longitudinal susceptibilities in zero 
field diverge. This is the inverse critical temperature of the system. For /5 > /5c, we find 
that the condition > is not satisfied on the whole interval < m < 1. Instead, a 
certain temperature-dependent value ms{(3) exists, such that yj_ vanishes for m = irts and 



becomes negative for m < ms- Since from Eqs. (|^), (|27|) , (pq) y± = implies H = 0, ms 
is just the spontaneous magnetization at the temperature considered. For m < nis not 
only does the transverse susceptibility behave unphysically, but the PDE itself becomes 
intrinsically unstable and cannot therefore be integrated numerically. Hence, we excluded 
the interval (0, ms) from the domain of integration. Below the critical temperature the 



boundary condition (Q) is then replaced by 

y^{m = ms,P)=Q. (46) 

The spontaneous magnetization ms is determined by checking the sign of y± as the 
integration with respect to /5 goes on. As soon as y_[_ becomes negative at a certain 
m = mo, the spontaneous magnetization curve is widened by setting ms = mo + Am, 
Am being the step of the magnetization grid. 

The numerical integration of Eq. (|42|) has been carried out with an inverse temperature 
step A(3 initially set at lO-'^-lQ-^ This can be further reduced as the critical region is 
approached so as to get very close to the critical temperature, and then gradually expanded 
back to its initial value. The number of mesh points in the magnetization grid ranged 
from ~ 10^ to ~ 10^. 



3 Results for the simple version 

In order to assess the reliability of the approximate closure (|3^), we considered first the 
case of a linear chain of three-dimensional spins, which can be solved exactly for zero 



magnetic field The SCOZA predictions agree very well with the exact solution, as 
can be seen from Figs. 1 and 2, where the results for the susceptibility and the specific 
heat are shown. We recall that this one-dimensional model has no phase transition at 
non-zero temperature. 

We then turn to the three-dimensional lattices. In this case the SCOZA does predict 
a phase transition at a certain temperature, which can be determined by locating the di- 
vergence of the isothermal susceptibility in zero field. The critical temperatures predicted 
by the SCOZA for the simple cubic (SC), body-centered cubic (BCC), and face-centered 
cubic (FCC) lattices are reported in Table I, together with the results of extrapolation of 
high-temperature expansions [|1^, ^ and finite-size scaling (FSS) analysis performed on 
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Monte Carlo simulations |^. The agreement is very good; the error is at most 0.55%. 
The critical temperature predicted for the SC lattice by the renormalizat ion-group based 
quantum hierarchical reference theory (QHRT) in the classical limit p2[ is also shown for 
comparison. Furthermore, the Table shows the critical internal energy per spin Uc- The 
SCOZA value of Uc is actually the same for any theory where the OZ approximation (|19D 
and the normalization condition (|T^ are satisfied. From Eq. ( pH]) with m = H = 0, 
z\\ = z± = 1 one has Mc,scoza = —<l[P{^) — l]/[2-P(l)]- We note that this result depends 
only on the type of lattice structure, so that for a certain lattice the SCOZA predicts 
the same Uc for any spin dimensionality. The agreement of the SCOZA result with the 
data from high-temperature series expansions appears to be only modest, and sensibly 
worse than in the Ising or D = 1 case |1|]. This is somewhat surprising, since the SCOZA 
is known to be exact in the D — > oo or spherical-model limit 0, so one might have 
expected very good agreement for intermediate D as well. Since the SCOZA prediction 
for Uc is not affected by the approximate closure relation (0), such a behavior seems 
to indicate that as the spin dimensionality is increased, the OZ approximation does not 
become more and more accurate in a monotonic fashion, but a "lowest accuracy" for some 
intermediate D might be reached. In Fig. 3 the reduced isothermal susceptibility in zero 
field above the critical temperature x\\ = X± and the correlation length ^ (BCC lattice) 
are compared with the results from extrapolation of high-temperature expansions . In 
both cases the agreement appears to be very good, even in the near-critical regime when 
^ ~ 10^. As stated in the Introduction, our numerical procedure for the D-vector model 
SCOZA gives for the critical exponents the spherical model values 7 = 2, z/ = 1, /3 = 1/2, 
a = —1, rj = 0, where the usual notation for the critical exponents has been used. In 
Fig. 4, the susceptibility in the critical region as a function of the reduced temperature 
t = (T — Tc)/Tc has been plotted on a log-log scale. We have also reported the effective 
critical exponent 7efr = —d{log x) / d{logt) , which shows that for the Heisenberg model the 
onset of the asymptotic critical regime as given by the SCOZA takes place for t ~ 10~^. 
This is in sharp contrast to the Ising case, also shown in the Figure, whose 7eff is affected 
by a very strong crossover, remaining very accurate over a quite wide range of t and sat- 
urating to its inaccurate asymptotic value only for t ^ 10^^. The specific heat ch in zero 
field has been plotted in Fig. 5 together with a closed-form approximant 0]. Because of 
its mean-spherical critical behavior, the SCOZA predicts that at the critical point ch will 
saturate at a finite value, which for the FCC lattice to which the Figure refers is given by 
Cpeak/^B = 3.03. For the Heisenberg model ch is indeed finite at the critical point. The 
approximant reported in the figure gives Cpeak/^B = 10, but one has to take into account 
that this quantity depends very sensitively on the value of the critical exponent a used 
in the approximant, which in Ref. was taken equal to a = —1/16. As reported in the 
same reference, using a = —1/8 shifts the peak to Cpeak = 5.7, which is probably a better 
approximation of the correct value as this choice of a is much closer to the best currently 
available estimate a ^ —0.122 [^. The SCOZA predictions for the magnetic structure 
factor r(k) in zero field (BCC lattice) are shown in Fig. 6 at two different temperatures, 
corresponding to reduced temperatures t ~ 2.2 ■ 10~^ and t ~ 1.4 ■ 10~^. The agreement 
with the approximant M is very good. 
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Let us now consider the results below the critical temperature. The spontaneous 
magnetization curve of the SC lattice is reported in Fig. 7, together with the results 
of the QHRT |2| and MC simulations [|]. In Fig. 8 the curve of the FCC lattice is 



compared to a Fade approximant. We see that in both cases the agreement is relatively 
poor, with the SCOZA results looking rather mean-field like. Actually, in the temperature 
range shown in the Figures, the main differences between the SCOZA and the mean- 
field phase boundaries can be accounted for by the shift in the critical temperature, so 
that they look very close to each other, provided they are plotted as a function of their 
respective rescaled temperatures T/Tc. The effective exponent P^s given by the SCOZA 
is plotted in Fig. 9. The critical behavior of the SCOZA phase boundary for the spin 
systems with continuous symmetry considered here is different from that found in the 
Ising or simple fluid models. In that case, the SCOZA critical exponents below do 
not coincide with the spherical ones found above T^, and the exponent /3 is given by 
/3 = 7/20 = 0.35, which gives a very accurate coexistence curve. Here instead we find 
spherical exponents both above and below Tc. We have then jS = 1/2, as in the mean- 
field approximation, so that the spontaneous magnetization curve close to the critical 
temperature is considerably narrower than the correct one, for which f3 ~ 0.36 On 
the other hand, is seems unlikely that this mean-field like behavior alone can account for 
the lack of global accuracy of the SCOZA phase boundary. In fact the QHRT yields for 
f] the non-classical value f3 = 0.41, but nevertheless is no more accurate than SCOZA as 
far as overall agreement with the MC coexistence curve is concerned. An improved closure 
of the SCOZA equation expected to deliver more accurate results for the phase diagram 
is discussed in the following Section. Here we observe that, despite its modest overall 
accuracy, the SCOZA phase boundary becomes exact in the low-temperature limit. In this 
regime one can use spin-wave theory, which gives for the spontaneous magnetization 



ms ~ 1 - -P(l) k^T + 0{T^) . (47) 

The spin-wave result for the spontaneous magnetization is compared to the SCOZA and 
to the mean-field predictions in Fig. 10. We see that the SCOZA predicts both the 
linear dependence on T for T — and the correct slope of the 0{T) term. Another 
point which deserves some attention concerns the behavior of the susceptibility across 
the spontaneous magnetization curve. For spin systems with continuous symmetry, we 
expect the longitudinal susceptibility to diverge not only at the critical point, but on the 
whole coexistence curve. Within the SCOZA, the coexistence curve should then merge 
with the spinodal. This is indeed what happens, as shown in Fig. 11, which reports the 
SCOZA inverse longitudinal susceptibility 1/X|| along two different sub-critical isotherms. 
In Fig. 12 several super- and sub- critical isotherms are compared to MC simulation 
results 0]. The discrepancy between the SCOZA and the MC data already observed in 
the coexistence curve affects the sub-critical isotherms for H ^ 0, especially for T near 
Tc. The agreement rapidly improves as H increases. 

In order to assess to which extent the results discussed here depend on the ansatz (P^), 
we also solved numerically the SCOZA PDF with a different closure relation, namely that 
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of Eq. (^21). We find that the thermodynamics and the phase boundary are hardly dis- 
tinguishable from those already shown, the critical temperature actually getting sligthly 
worse: for instance, for the SC lattice one has Pc = 0.6978, while the previous approxi- 
mation (see Table I) gave pc = 0.6968, and the "exact" resuh is §, |HTT[ /J^ = 0.6930. 
Moreover, our algorithm failed to find solutions below a certain sub-critical temperature. 
We do not know whether this is a purely numerical problem, or instead it is the sign of a 
deeper disease which may affect the closure (|32D. In this respect, we recall that in such an 
approximation (see Sec. II) the mean-square fluctuations of the spin in the transverse and 
in the longitudinal directions do not behave correctly on the coexistence curve, while the 
correct behavior for m = 1 is enforced by the boundary condition (^) . It is then possible 
that at low temperatures, as the spontaneous magnetization approaches its limiting value 
ms = 1, the change in behavior across the domain becomes so abrupt as to prevent the 
existence of a solution. Since the above closure does not introduce any improvement with 
respect to the approximation (0) that we have used so far, we have not investigated this 
point any further. 

4 Inclusion of induced spin- wave correlations 
4.1 New closure 

We now seek to develop an improved version of the SCOZA, guided by the results for 
the simple version that we have just summarized. The most serious deflciency of that 
version is its spherical-model like coexistence curve, which considerably differs from what 
is expected to be the correct result. Thus we have to look for a possible leading correction 
that might remedy this. What we intend to do is to explicitly take into account the 
effect of the transverse correlations on the longitudinal susceptibility. As we will see 
below, this will naturally give rise to a diverging longitudinal susceptibility x\\ the 
coexistence curve by enforcing y\\ = whenever y± = 0. On the other hand, if the effect 
of transverse correlations is not properly accounted for, a diverging x\\ at coexistence can 
be obtained only by heavily distorting the longitudinal part of the correlation function, 
with a corresponding infuence upon the internal energy and thus the SCOZA results. 
This is most probably what happens in the simple version that we have considered so far. 
For magnetization m — 1 one approaches Gaussian-model behavior and these induced 
correlations can thus be obtained exactly by a spin-wave calculation. For general m, we 
can turn to the analysis based on 7-ordering in order to obtain their form, where 7 is the 
inverse range of interaction. 

The induced correlations become long-ranged near the critical point and at phase 
coexistence. Thus their major influence will be upon the longitudinal susceptibility x\\- 
Representing their contribution integrated over space by A we can write, for the longitu- 
dinal susceptibility. 
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where the new quantity yi represents the contribution of the MSA-hke part of the corre- 
lation function. Solving with respect to yi one then has 



„ = ^i^. (49) 

and there will be a quantity zi related to yi by Eq. (^). The direct contribution these 
induced correlations will make to the internal energy will be minor and less crucial ac- 
cording to 7-ordering, so we neglect it here for simplicity. The quantity yi is then used 
instead of y\\ in Eq. (|38| ) to evaluate the longitudinal part of the internal energy. This 
requires having an expression for A. The analysis is most conveniently given in the dia- 



grammatic language developed ||2^ to describe the 7-ordering used by H0ye and Stell. In 
that language, the correlation functions are built up of chain bonds (representing chains 
of pair interactions) and hypervertices (representing the short-range components of the 
correlation functions). The leading contribution beyond the MSA to the pair correlation 
function r(r) of a fluid system has the diagrammatic form 



r(r-r ) - r (r-r ) 




where the triangles represent the three-body correlation function of the reference system 
with purely short-range interactions (including self-correlations where two or all three of 
the particles can be identical), the white circles (or "root points") represent the coordi- 
nates ri, r2, and the full circles represent coordinates that are integrated over. The zigzag 
chain bond C(r) (not to be confused with the direct correlation function c(r)) is defined 
by 

C(rj-r2) = cAAA/Wb = o - - - <3 + o - - - -• • -o + o -• • -• » <3 + . . . . 



where the solid and dashed lines represent respectively the two-body correlation function 
of the reference system rR(r) and the perturbing function —(3ip{Y). In Fourier space one 
then has 

(7(k) = -^^^^j . (50) 

i + /3^(k)rR(k) 

The above result can be generalized to spin systems by utilizing Onsager's recipe for 
regarding different orientations as different species of a mixture l^^j. In particular, the 



reference system consists of non-interacting spins, and is mapped into a multicomponent, 
but simple, hard-core lattice gas whose n-body correlation functions T^^iri . . .Yn) are 
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non-vanishing only when the coordinates ri . . . r„ coincide. In particular, the three-body 
hypervertex has the form 

rfjVbr2,r3) = b^nra^rara • (51) 

(3) 

In the multicomponent case, the Fj^ as well as the chain bond should be summed over 
the various species. One finds 

where fig, the chemical potential has the form sH + const. Thus the sum on the left- 
hand side can be replaced by a differentiation with respect to H, which implies that the 
resulting rj^"* can again be represented by the one-component form Eq. (51). For the 
multicomponent system considered here, b is a three-index tensor with components bijk, 
where each index can assume D values corresponding to the longitudinal and each of the 
D — 1 = d transverse components of the spin, b^k is determined by the sum rule 



^ijk = ^ L Tr^' (ri,r2,r3) 



(53) 

H«=H,H, =0 



where is the number of lattice sites and c<jr is the grand potential per site of the 
reference system. This is the analog of the compressibility rule (|T^), (T^) for the two- 



body correlation function. We recall that for the reference system the latter gives x\\ = R\\ 
Xl = R±, where the quantities R\\, R± have been defined in Eqs. ([I5|), (|16D. According to 
Eq. (|12|), the contribution A in Eq. (^) is obtained by integrating over space the leading 
contribution to r(r) — rMSA(r) shown above, where each of the white circles corresponds 
to a longitudinal root point, and the full circles involve both a spatial integration and 
a summation over the component index. We observe that the matrix of the two-body 
correlation function (see Eqs. (^, (|T3)) and that of the nearest-neighbor spin-spin 
interaction s^sl, are diagonal. Eq. (^) then implies that this is true also for the matrix 
of the chain bond Cjj(r), i.e., the zigzag bond is non- vanishing only if its endpoints 
correspond to the same spin component. It is then readily seen that in performing the 
summation over the indexes, the only two terms which contribute to A are those in which 
both the three-body hypervertices are either of the type or of the type By 
setting ^ = b and applying Eq. (^31) we find 



R\\R'^ = —{R\\-Ri_), (54) 



where the prime denotes differentiation with respect to m, and Eq. (|^) has been used. 
We note that, since R\\ — R± = 0{m?) for m — 0, 6 vanishes for m = 0. In the limit 
m — > 1, = -R||-Rj| can be disregarded with respect to b as then one has 6 ~ (1 — m)^, 
6|lllll ~ (1 — mY (see the following section). We have then 

^ = ]-db'^B, (55) 
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where the factor 1/2 is the graph symmetry factor with respect to the interchange of the 
two chain bonds, d is the number of the transverse components, and B is given by 



(56) 



C±{v) being a short-hand notation for C±^{t). The chain bond follows from Eq. (|50|), or 
equivalently from the MSA-like form of the correlation function, whose Fourier transform 
is given by Eqs. (p 



(P0|). The transverse part is then 
f.(k) ^' 



1 - Ri. (4 + gc^7(k)) 



(57) 



where Cq can be regarded as representing the perturbing interaction inside the hard core. 
With = /3 as in the MSA, i?_L is clearly the two-body hypervertex of the reference 
system. With ^ j3, this interpretation can be kept as the difference — j3 is small 
and will become even smaller when including the new contribution A. We then write 



Along with Eq. 



f±(k) = i?x + i?lC'±(k), 
this yields for the chain bond in Fourier space 



C'x(k) 



7(k) 



1 - ^±7(k) 



{51 



(59) 



Note that since T±{r = 0) = R± we have C±{r = 0) = 0, which also follows from Eq. ( p9D 
using the definition (p3D for the integral P{z). Now the integral of the double chain bond 
can be evaluated in a straightforward way by Fourier methods: 



B 



(2^ 



1 



RlP{zi.) 



where we have set P'{z) = dP{z)/dz. Inserting Eqs. ( ^41) and 
the sought quantity A to be used in Eq. 



[^±P'(^±)-P(^±)(P(^±)-1)] , (60) 
I) into Eq. (|55D yields 



4.2 Spin-wave limit 

As m — i> 1 the transverse modes represent harmonic oscillators that can be solved exactly 
by Gaussian integrals. With the added contribution A of Eq. ( ^8]) the SCOZA will 
recover the exact spin-wave solution also for the longitudinal susceptibility as m — > 1. I.e. 
thermodynamic self-consistency is immediately fulfilled with cl = Ci = f] in this limiting 
case also at phase equilibrium where P'{z±) oo. 

To show that SCOZA has this behavior, we first evaluate R± as m — 1. By setting 
u = d/2 we have 

m = (4) = {^l-d{s^y) ~ 1 - u{{s^r) = l-uR^ (61) 
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so that in this hmit 

flx = i^. (62) 

If we set Ci = P, Eq. ( pS]) then gives 

q (3{1 - m) = uz±P{z±) . (63) 

We also note that in the same hmit R\\ = {{sj — (sl!))^) = (1 — m)^/z/, so that it can be 
neglected compared with R±. Eqs. (p4[) and (p^) then imply 



6^-<i^. (64) 



Z/2 



This together with Eqs. (^), (0)' O) yields 



A = (1 - "^)'^;p^ [z±P\z^) - P{z^){P{z^) - 1)] . (65) 



In order to obtain the longitudinal susceptibility from Eq. ( ^8]) one also needs knowledge 
of?//. Eq. ( p5D with cf = /? and z = zi implies zi — > qf3R\\ q(3{l—mY/h' (i.e. P{zi) 1), 
so that Eq. (|36D gives 

= 1 - ^ 1 . (66) 
Substituting Eqs. (|65D and (|6^) into Eq. (|48|) we obtain 



1 ,2 ^±P^(^±)+P(^±) 

XII = -(1--) • (67) 

This can be compared with the result obtained by differentiating the transverse suscep- 
tibility with respect to m. For m 1 only the transverse correlations contribute to the 
internal energy, and the SCOZA equation (^) reduces to 

1 A (±] = -1 _ A („R^^^M^] = _i _ J^p'^,^)^ (68) 

q dp \X± / 9m \ z±P{z±) J qjS dm 



where again Eq. (|25|) with cj- = (3 has been used. In the same limit we then also have 
with Eqs. (0) and (|§ 

y± ^ Piz^) - q(3R± (69) 



so that from Eqs. (|27|) and (^ 



The partial derivatives of z± which appear in Eqs. (|68D, (|TOD are obtained by differentiating 
Eq. (|63|) with respect to /5 and m. This gives 

dz± ^ q{l-m) 

d(3 i^[z^P'iz^) + Piz^)]' ^ ' 

^ = ^1 . (72) 
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By substituting Eq. (ffll ) into Eq. ( [70D it is readily seen that the consistency condition 
Eq. is satisfied. Finally by using Eqs. (§, d^^), (|T2|) the longitudinal 

susceptibility becomes 



in agreement with Eq. (^) . Thus the contribution A in Eq. (^Sf) fully accounts for the spin 
waves in the limit m — > 1, and SCOZA yields full thermodynamic self-consistency with 
(\ = c\ = (5 m. this limit. Without this term Eq. ( |67D reduces to X|| = (1 — mY/iJ = R\\ as 
a consequence of the assumption c\ = (3 that leads to Eq. (|65|) . However, SCOZA, that 
does not "know" about this assumption, still requires Eq. ( |73| ) to be fulfilled (with c]^ ~ /3) 
on the basis of the MSA-type expression (|18|), which implies X|| = ~ R\\{cl + ?cf)]- 

At phase equilibrium where X|| oo since P'{z±) — > oo this implies that Cg + cf — oo 
as -R|| ^0, so that c\ — (3 is no longer small and the longitudinal part of the correlation 
function will be heavily distorted. 

5 Conclusions 

In this paper a simple version of our SCOZA is studied, in which the parallel and transverse 
components of the onsite spin-spin correlation R are assumed to be those of the high- 
temperature susceptibility. This is found to be quantitatively accurate in its assessment 
of Tc but not in several other respects. This treatment neglects the effects of the transverse 
correlations on the longitudinal susceptibility. We then develop a version of SCOZA in 
which the neglected correlations are included to lowest order in their 7-ordered expansion. 
This guarantees that spin waves are well-accounted for in a way that becomes exact in 
the m ^ 1 limit. The quantitative assessment of the resulting version of SCOZA will 
require a significant computational effort and is not attempted here. We defer it to our 
next paper on the Z)-dimensional spin model. 
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TABLES 



lattice 


/^cSCOZA 


/5c,QHRT 




■^c.SCOZA 




SC 


0.6968 


0.7047 


0.6930*^ 


-1.0216 


-0.987 ±6^ 


EGG 


0.4887 




0.4868" 


-1.1289 


-1.092 ±8^ 


FGG 


0.3160 




0.3148^ 


-1.5379 


-1.4808 ±30^ 



Gritical parameters of the classical Heisenberg model on a SG, BGG, and FGG lattice. 
The inverse critical temperature (3c and the critical internal energy per spin Uc predicted 
by the SGOZA are compared with the "exact" /3c,ex, ""cex determined by extrapolation 
of high-temperature expansions and by MG simulations. Energies and temperatures are 
in units of the interaction strength J and of J/k-Q respectively, a: results from high- 

, ^. h\ from high-temperature 
For the SG lattice, 



temperature expansions |TT| and MG-FSS simulations 
c: from high-temperature expansions I 



expansions 



10 



the inverse 

critical temperature predicted by the quantum hierarchical reference theory (QHRT) 
in the classical limit is also reported. 
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Figure 1: Reduced isothermal susceptibility x iii zero field of the classical Heisenberg 
model on a one-dimensional chain as a function of the dimensionless temperature k-^T/ J, 
J being the spin-spin coupling strength. Solid line: SCOZA. Crosses: exact result pl| . 
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Figure 2: Specific heat per spin ch in zero field of the classical Heisenberg model on a 
one-dimensional chain as a function of the dimensionless temperature k^T/J. Solid line: 
SCOZA. Crosses: exact result [^]. 
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Figure 3: Reduced isothermal susceptibility x and correlation length ^ (b) of the 
classical Heisenberg model in zero field above the critical temperature (BCC lattice). ^ is 
in units of the lattice spacing /. Crosses: SCOZA. Solid line: approximant M. 
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Figure 4: Reduced isothermal susceptibility x iii zero field (a) and effective critical 
exponent 7efr (b) of the classical Heisenberg model (BCC lattice) as a function of the 
reduced temperature (T — Tc)/Tc. Solid line: SCOZA. Dashed line: approximant 0. 
Dotted line in panel (b): SCOZA result for the Ising model 
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Figure 5: Specific heat per spin ch in zero field of tfie classical Heisenberg model above 
the critical temperature (FCC lattice). Crosses: SCOZA. Sohd line: approximant 0. 
The arrow marks the height of the peak of the SCOZA specific heat Cpeak — 3.03. 
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Figure 6: Magnetic structure factor r(k) in zero field of the classical Heisenberg model on 
a BCC lattice along the direction kx = ky =kz. The norm of the wavevector |k| is in units 
of the reciprocal of the lattice spacing /. The results shown are for inverse temperature 
(3 = 0.4 (a) and (3 = 0.48 (b). Crosses: SCOZA. Solid line: approximant §. 
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Figure 8: Spontaneous magnetization curve of the classical Heisenberg model on a FCC 
lattice. Solid line: SCOZA. Squares: Fade approximant 0. 
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Figure 9: SCOZA effective exponent fJ^s for the curvature of tlie spontaneous magneti- 
zation curve of the classical Heisenberg model (SC lattice) as a function of the reduced 
temperature (Tc — T)/Tc. 
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Figure 10: Low-temperature behavior of the spontaneous magnetization curve of the clas- 
sical Heisenberg model (SC lattice). Solid line: SCOZA. Dotted line: mean-field theory. 
Squares: MC simulation data 0]. Dashed line: spin-wave theory result to 0(T) [0. 
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m 

Figure 11: Reciprocal of the reduced longitudinal susceptibility l/x|| of the classical 
Heisenberg model (SC lattice) outside the phase boundary. The SCOZA results are shown 
for two different values of the inverse temperature /5. 
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Figure 12: Magnetization isotherms of the classical Heisenberg model (SC lattice). The 
magnetization m is plotted as a function of the magnetic field H for the values of the 
inverse temperature (3 shown in the figure. Lines: SCOZA. Squares: MC results [Q. 
Dashed and solid lines refer respectively to super- and sub-critical isotherms. 
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